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ABSTRACT 



We study the formation of the Calcium II infrared triplet lines 8498A, 8542A and 8662A, in the accreting magnetospheric flows of Classical 
T Tauri stars (CTTS), and present a grid of models for a large range of magnetospheric conditions. We apply our models to the interpretation 
of multi epoch observations of the CTTS DI Cep. We find that these lines form in the magnetospheric infall and that the variability of the Call 
triplet lines in DI Cep can be explained in the context of changes in the mass accretion rate/temperature of the accretion column gas. 
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1. Introduction 

T he concept of magn etically controlled accretion, introduced 
bv lGhosh et al. ( 1977 ) for neutron stars, was applied to T Tauri 
stars b y lUchida & Shibatal d 1985b . ICamenzindl dl990l) . iKonigl 
(1991)' and lShu et al.l (11994). From the interaction between the 
magnetic field lines and the disk, by transporting angular mo- 
mentum outwards, it is possible to get the braking torque nec- 
essary to spin down T Tauri stars (TTS). Magnetocentrifugally 
driven winds can also remove angular momentum from the disk 
as it recei ves stellar angul ar momentum via the magnetospheric 
coupling ( Shu et a ll 19941). T hese theories can explain observa- 
tions bv lBouvier et al.l d 1993b and lEdwards et alj (11993b . which 
show that Classical T Tauri Stars (CTTS), which accrete from 
disks, rotate slower than the non-accretors Weak T Tauri Stars 



(WTTS). However. iMatt & Pudritzl d2004h argue that the disk 
braking torque may not be enough to spin down T Tauri stars 
at all. 

In the con t ext of CTT S Irlartmann et alj d 1994b and 
Muzeroll e et al. ( 1998 , 200 lb developed axisymmetric models 
of spectral lines formed in the magnetospheric accretion flow 
towards the star. They showed that the observed characteris- 
tics of the line profiles, such as (1) slightly blueshifted emis- 
sion peaks, (2) blueward asymmetries, and (3) redshifted ab- 
sorption components/inverse P Cygni profiles, could be well 
described by the models. The blue/UV excess emission com- 
ing from the high velocity materi al that shocks onto the star 



was also successfully modeled by Calvet & Gullbringl (1998) 



and iGullbring et al. (2000). The magnetic field line intensi- 
ties necessary to disrupt the ma terial from the di s k are con- 
sistent with measurements from Johns-Krul l et al. d 1999b and 
Johns-Krull & Valentil d2000b . 



Most modeling of lines supposedly formed in the magne- 
tospheric accretion flow has been focused on hydrogen and 
sodium lines. But other tracers of the complex environment 
around TTS are also observed, such as the infrared triplet lines 
of Call 8498A, 8542A, 8662A, or the Hel lines 10830A (ac- 
tually a triplet) and 5876A, among others. We have started a 
modeling effort to simultaneously reproduce observations of 
these different lines, in order to restrict parameter space and 
to test the consistency of the model. These parameters include 
the mass accretion rate, the temperature of the gas and geom- 
etry. As a step in this direction we develop a study of the Call 
infrared triplet. 

The Call infrared triplet line profiles from CTTS can be 
classified into three distinct groups, according to their shapes: 
only narrow component profiles (NC), only broad component 
profiles (BC), profiles with simultaneous NC and BC. When 
only the NC is present, the profiles are very similar to those 
of WTTS, which display enhanced narrow emission at the 
cente r of the photospheric absorption line (Batalh a & Basri 



1993; Batalha et al. 1996; Muzerolle et al. 1998). BC emis- 



sion lines hav e profiles extending to velocities in excess of 
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100 Km s~' dHamann&Perssonl ll99l iBatalha et al] 1 19961: 



Muzeroll e et alj Il998b. The spectral samples presented i n 



Hamann & Perssonl d 1992b ; (sample 1) IBatalha et ail dl996b : 
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(sample 2) and iMuzerolle et al.l d 19981) ; (sample 3) give a full 
overview of these possibilities. 



2. Models 

2.1. Line profile models 

We can obtain further insight into the infrared triplet line char- 
acteristics by comparing them to the more commonly studied 
H a line. Using, for example, sample 1 and sample 2 one can 
compare the Call 8542A profile with that of H a for about 43 
CTTS. The main conclusion of this qualitative analysis is that 
only for very few cases H a and the Call line exhibit similar- 
ities both in shape and extension of the line wings. The Call 
line is in general less featured. Using other lines of the in- 
frared triplet would lead to the same conclusions, since there 
are no significant structural differences among them. Profiles 
with blueshifted emission, slight redshifted absorption, and su- 
perposed blueshifted absorption can appear in both H a and the 
Ca II infrared triplet lines. 

The Call lines show a p eculiar pattern, found in m any T 
Tauri and Herbig AeBe stars (Ha mann & Perssonlll99d) . If the 
8498A, 8542A, 8662A lines were optically thin, their strength 
would follow the 1:9:5 ratio according to their (gf) values of 
0.055,0.49,0.27. The observed ratios a re usually close (but 
not eq ual) to the optically thick limi t dHerbig & Soderblom 



not equal) to tne optically tnicK limi t OHerbig & bo aeroiom 
1980t Hamann & Perssonl 1 1989i Il990l) which corresponds to 



the ratio 1:1:1. The more extreme observed peak intensity pat 
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(Hamann & Persson 1990, 1992). There is a tendency of the 



8498A line to be the narrower line, specially when it becomes 
stronger, w hile the width of 8662A line is sim ilar to that of the 
8542A line dHamann & Perssonlll99dll992h . 

A chromospheric origin of the NC p resent in Call line 
profile s of CTTS has been pr opo sed bv lHamann & P erssonl 
dl992h . iBatalha & Basri (1 19931) and lBatalha et alJ dl996l) . with 
the enhancement relative to WTTS being driven by disk ac- 
cretion phenomena ac t ing indirectl y on the chromosphere 
dBatalha & Basrilll993l) . lBasril dl990l) argue that the BC could 
be explained by the sum of emission arising in turbulent 
regions of different velocities, which would also account 
for the narrower than gaussian c entral peaks seen in many 
CTTS. lHamann & Perssonl d 1992b also suggested that the BC 
Call lines formed by turbulence in an extended envelope. 
Gullbring et al. d 19961) argued that the BC could be formed 
in the high density post shock region of the accretion flow. 
According to lCalvet & Gullbring! ( 19981) . the temperature range 
immediately after the shock is around lO 5 ^ - 10 K which 
is probably too high to leave enough Call available. The ve- 
locity is typically of the order of 100 Km s _1 in this region 
which is less than the velocity extent of the lines commonly 
observed. Deeper in the postshock where the temperatures are 
lower, the velocity has decreased even more. On the other hand, 
Muzerolle et al.l dl998l) showed that models based on magneto- 
spheric accretion infall could reproduce the observed profile of 
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Fig. 1. Call Grotrian diagram showing the five energy levels 
plus continuum used in our models. The radiative transitions 
considered are represented by the dashed lines. 

the BC observed in BP Tau, and were consistent with the pro- 
files of other species. 

In this paper we explore Call triplet line models further, 
by spanning a larger range of magnetospheric physical condi- 
tions. Although there is evidence for non-axisymmetric accre- 
tion a s a way to explain the observed variability and the hot 



spots dBertout et all 1 9881; iBouvier et all 19931: IMuzerolle et al 
2001; iRomanova et al.l I2.QQ41) , we use an axisymmetric model 
for simplicity. We com pare the model result s with observations 
of DI Cep collected by lGameiro et al.l d2005l) . hereafter GFP05. 
We discuss which values of the parameters yield better fits, 
trying to constrain them and comparing to those derived from 
other work. We will also attempt to explain the observed line 
variability of DI Cep by changing the crucial magnetospheric 
accretion parameters,!. e. the mass accretion rate and tempera- 
ture of the gas. 

We developed models of the Call infrared triplet lines 
8498A, 8542A and 8662A. The model atom consists of a five 
level ion plus continuum, with levels As 2 S 1/2, 3d 2 .03/2,5/2 and 
Ap 2 P\/2,3/2 (se e figureED- The l i ne pro files are computed based 



on the work of lHartmann et al.l d!994l) . We present in this sec- 
tion a brief review of the method. 

In the adopted model, the disk is truncated by the mag- 
netic field, which has an assumed axisymmetric dipolar geom- 
etry and funnels the material towards the surface of the star 
where it arrives with near free fall velocity, merging onto the 
photosphere through a shock, si mulated by a ho t ring. The 
star is not rotating: as shown by IMuzerolle et alJ (2001) for 
the majority of the stars with slow rotational velocities, the 
effect is negligible in the final profiles. The numerical code 
CV, developed by lHartmann et al. ( 1994) to solve the radia- 
tive transfer in this kind of system, is used to calculate the 
source functions and radiation fields for the calcium lines. This 
code is based on the exte nded Sobolev approximation from 
Rvbicki & Hummerl (1 1978b which allows for a generic veloc- 
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ity flow geometry in relation to the standard Sobolev approxi- 
mation. We then solve the statistical equilibrium equations to 
find the atomic level populations. The atomic data was col- 
lected from different sources: the collisional cross sections 
from PANDORA program; the A^fa, energy level values were 
taken from the online Atomic Spectra Database V3.0 available 
at the US National Institute of Standards a nd Technolog\o 
and th e photoionization cross-sections from IShine & Linskv 
( 11974 . 

Initially we set the input parameters of the model, i.e., the 
mass accretion rate (M acc ), the geometry of the magnetosphere 
specified by the range of radii (ro - rf) which corresponds to the 
range where material is loaded from the disk into the magneto- 
sphere, the temperature structure of the accretion flow charac- 
terized by its maximum value (T gas ) and the inclination (i), The 
effective temperature of the star (T e ff) has also to be specified 
as well as the mass (M) and radius (/?). The temperature of the 
shock (T r ) can be derived from the previous parameters in light 
of the magnetospheric equations but we will assume it to be a 
8000K blackbody, unless otherwise indicated. 

With these parameters set we use CV to generate the hydro- 
gen and electron density structure across the magnetosphere. 
The temperature structure given the maximum v alue (T gas ) is 
calculated as described in lHartmann et al.l (1 1 994 .Considering 
the ionizing radiation field at a given point of the magneto- 
sphere as the sum of a Planck contribution from the star and 
the ring multiplied by an adequate dilution factor, 



J v — W starBy(T star) W r i n gB y (T r i n n) 



(1) 



with the Jj mean intensities of the 5 level system allowed tran- 
sitions to be equal to the local Planck function, we can make 
an initia l gues s for the populations in each level, using the 
Mihalai d 19781) formalism, and a standard LU decomposition 
routine. We then use these populations to retrieve new J, from 
CV, which we include in the statistical equilibrium equations, 
with the purpose of calculating new populations, repeating this 
procedure until the populations converge. Finally, the line pro- 
file is calculated summing all the emission from the re sonant 
surfaces at a given frequency (see Hartmann et alj[l994 ). 

In this section we explore the results from different input 
model parameters and initially describe some basic quantitative 
properties. 

2.1.1. Density 

In Fig. [2] we show the densities along the magnetosphere for 
a given set of parameters. The assumed abundance of Ca rel- 
ative toH is 2.37xl0~ 6 (similar to the solar system value in 
Cox 2000). From setting the geometry and M acc alone, we get 
in (a) a density behavior along the magnetosphere regulated 
by mass conservation. As the material is free falling, the in- 
crease inwards is due to the confinement of the material in a 
smaller area. As the velocity decreases rapidly near the disk 
and the area of the flow is more or less the same, the density 
rises. Comparing the density ratio of Call/Calll in LTE and in 
NLTE, shown in figure|2](b), we can infer that the Call density 



http://physics.nist.gov/PhysRefData/ASD/index.html 



is much higher in the models than expected in LTE and the be- 
havior is also different. The only continuum sources considered 
are the star and the accretion shock which have typical tem- 
peratures lower than those in the magnetosphere. Radiative re- 
combination is therefore more important in the magnetosphere 
than photoionization in rel ation to a LTE case (defined as in 
equations 5-66 and 5-67 of lMihalasI 1978I) . so the Call density 
increases by some orders of magnitude. Were the ratio between 
radiative recombination and photoionization constant through- 
out the magnetosphere and the variation of the ionization state 
with distance would be similar to LTE. The LTE behavior is 
just a consequence of the temperature profile (see Fig.[3]Planck 
functions), and to a much less extent the electron density, so 
that is why Calll is more important in the middle of the mag- 
netosphere. In our NLTE case the ratio of Call/Calll is flat- 
ter than in LTE. This happens because the dilution factors are 
high and rise fast towards the star, and so photoionization rises 
fast too. Far from the star, about > 1.9/?+ for the streamline 
that intersects the disk at 2.2/?+, the dilution factors decrease 
slowly and the gas temperature characteristic of radiative re- 
combination decrease fast, so the rise in the ratio Call/Calll is 
not strong. Near the base of the magnetosphere, near the disk, 
the decrease of temperature is so fast that the ratio decreases, 
towards the LTE values. Calll is the dominant state at g 1.8/?* 
and Call the dominant state afterwards (panel c). However, the 
Calll density rises near the disk, since photoionization and ra- 
diative recombination are more balanced again. Note that in a 
larger geometry than the one we show, a narrower flux tube 
near the star leads to higher total densities and electron den- 
sities. As a result, in panel (c) we would see the Call density 
increase towards the star, following a behavior similar to that 
of H density in Fig. [2] (a), but more attenuated because of the 
photoionization rise. Panel (d) shows the departure coefficients 
from LTE for the five levels. The high values of the departure 
coefficients result from the increased Call density in relation to 
LTE. The departure is higher for the lower levels and higher in 
the inner regions of the magnetosphere, which is evident in the 
continuum balance result shown in panel (b). 



2.1.2. Source Functions 

Figure [3] shows the behavior of the Planck and source function 
for the 8662A line along one streamline, for different param- 
eters. In (a) we use a smaller geometry in two different mass 
accretion rate regimes and in (c) a larger geometry, where we 
also study the effects of different gas temperatures. The ratio of 
the source to the Planck functions are plotted in (b) and (d). In 
the low M acc case in (a), we see a rise of the source function 
near the star. This happens because the collisional excitations 
from the metastable level 3d 2 Dj/2 to the level 4p 2 Pi/2 increase 
significantly due to the increase in the electron density. As we 
move further out, the source function decreases while the ef- 
fect of the electron density in the collisional excitation is dom- 
inant. At a certain point the source function rises again because 
the Call densities become rapidly higher and so the opacity 
of the H and 8662A lines increases. The escape probability de- 
creases significantly and the upper level remains overpopulated 
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Fig. 2. Densities across the magnetosphere. (a) Total densities of H (solid line) and Ca (dashed line) across five streamlines, that 
intersect the disk at 2.2, 2.4, 2.6, 2.8 and 3/?*. Ne is the electron density (dotted-dashed line), (b) Ratio between the density of 
Call and Calll in our NLTE models (dashed line) and in LTE (solid line), (c) Densities of Call (solid line) and Calll (dashed 
line), (d) Departure coefficients from LTE for the five levels in study. The departure coefficients are plotted for the streamline 



that intersects the disk at 2.2/?* 
-8.0, Geometry — 2.2 — 3R*, i 



R* = l.SR@, 
--60°. 



Mi, = 1.42M , T* = 5400/T, T r = 8000/C T gas = 12000/T, log{M acc ) = 



in relation to the lower one. Near the end of the magnetosphere 
the electron density drops very fast (much lower temperatures), 
and so the strong decrease in the collisional excitation makes 
the source function decrease. The temperature drop effect is 
therefore not present in the ratio of the source function to the 
planck function plotted in panel (b). 

In the high M acc situation the increase in the Call density 
makes escape in the line more difficult, and the increase in 
the electron density allows more collisional excitations. These 
two effects produce a source function that is more thermal- 
ized. Using panel (c), we can compare the high temperature 
case with the low mass accretion one of panel (a). Both cases 
only have different geometries. The rise of the source function 
towards the star is much more evident in the large geometry. 
This happens because the Call density rises towards the star 
in the large geometry (see section |2. Lit . As a result we have 
more radiative trapping in relation to the smaller geometry case 
which makes the increased collisional excitation towards the 
star more evident. The source function for the lower gas tem- 



perature regime does not have a different behavior. As shown 
in panel (d), it is less thermalized than in the high temperature 
case, because of the lower electron densities, since hydrogen is 
less ionized. 

In Fig. [4] we plot the source functions for the five lines 
arising in the Call system. We use the same parameters as in 
the low gas temperature case of Fig.[3](c) with the exception 
of a narrower geometry: 4.8-5/?*. Note that the 8662A line 
is more thermalized than in Fig. [3] (c), because narrower ge- 
ometries imply higher densities. One can see that the source 
functions for the five lines arising in the Call system are close 
to each other (valid to other sets of parameters). Near the star 
and close to the disk, the H and K source functions approach 
each other as well as the IR triplet lines. In other regions the 
source functions diverge, because the lines are less thermal- 
ized in the inner regions of the magnetosphere. The 8498A 
source function stays close to the 8542A source function since 
the lines share a commo n upper level. We also note the rela - 
tionS 42 /S 4 i ~ S52/S51 dShine & Linskvlll974l IZirkerll 1965b . 
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Fig. 3. The Planck (B) and source (S) functions behavior for the 8662A line along a streamline, for different parameters. In (a) 
and (b) we use the geometry specified by the radii 2.2-3/?*. B,S and S/B are plotted for the streamline that intersects the disk 
at 2.2/?*. A high M acc = 2 x l0^M Q yr~ l case and a low M acc = 1 X \Q^M yr' 1 case are shown. In (c) and (d) we used the 
geometry specified by the radii 4.2-5/?+. B,S and S/B are plotted for the streamline that intersects the disk at 4.2/?+, in a high 
T gas = 12000/^ and a lower T gas - 100007^ temperature regime. Other parameters as in Fig. [2] 



which means that the changes in the H and K lines source func- 
tions affect similarly the 8662A and 8498A line source func- 
tions. 



2.1 .3. Line profiles 

2.2. Exploration of model parameters 

In this section, we show line profiles for a range of mass ac- 
cretion rates, temperatures and geom etries. The behaviors w e 
find are similar to tho se indicated in lHartmann et alj dl994l) ; 
Muzerolle et al.l(l200ll) . although we do not consider the effects 
resultant from including line broadening mechanisms and con- 
tinuum opacity. We neglect radiative, Van der Waals and Stark 
broadening. Without a proper ray-by-ray treatment it is diffi- 
cult make a simple argument. However, the broadening will 
be much less important than for example in H a or Hp. For 
the exploratory nature of this paper therefore the assumption 
is adequate. As shown in Fig. [5] fixing the temperature, for de- 
creasing mass accretion rates, we find less emission due to the 
decrease of densities. Strong redshifted absorption appears be- 
cause of the growing importance of the background continuum 



relatively to the line source function (Mu zerolle et al.ll200ll) . 
When fixing the mass accretion rate and increasing the temper- 
ature the lines become more thermalized, resulting in higher 
line emission, and in more symmetric profiles. For the low mass 
accretion rates and low temperatures the model profiles do not 
have observational correspondence. For the highest mass accre- 
tion rates and temperatures, t he line does not go into absorption 
as in iMuzerolle et al. (2001) because of no continuum opacity 
treatment in our work. 

In Fig.[6]we show profiles for various geometries and val- 
ues of Macc- For better illustration, we use a high T gas for low 
M acc . We see more emission in the larger geometries than in 
the smaller ones, specially near the line center. The profiles are 
also narrower and more symmetric. In a large geometry, for the 
same inclination angle (60°), the line of sight crosses much less 
higher velocity material resulting in narrower profiles. Large 
geometries have higher emitting volumes, and due to the more 
extended flow and less curvature, the line of sight for this in- 
clination crosses more surfaces where the velocity projected in 
the observer's direction is nearly zero, resulting in increased 
emission, mainly near the line center. The increased symmetry 
results from less occultation effects from the star and the disk 
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Fig. 5. Normalized 8662A line profiles for different mass accretion rates and high temperatures of the gas. Fixed parameters: 
R+ = l.9R@, M-,, = O.8M , = 4000/T, T r = 8000/T, Geometry = 2.2 - 3fl*, i = 60°. 



(Muz erolle et al.l 120011) . The difference in intensity when we 
compare wide and narrow geometries resides in the fact that in 
the latter case we have higher densities across each line of sight, 
since the accretion rate is the same and we confine the material 
in a smaller volume. In some of the profiles there seems to be 
a NC on top of a BC however they are wi der than the typical 
NC widths found in C TTS of < 50Km s' 1 dBatalha et alJl996t 
Beristain et alhoOlh . 



In Fig. [7] we explore the inclination dependence. For high 
inclinations the profile is broader since our line of sight crosses 
more high velocity material. For low inclinations, a deep ab- 
sorption is seen near zero velocity that comes from the fact 
that our line of sight is parallel to the flow near the disk 
which is assumed to have a startup velocity of 10 km s _1 . 
Th e behavior is the same as in the hydrogen profiles modeled 
by lHartmann et al.l d 19941) . remarking that th e broadening ef- 



fects can smooth slightly the profile structure (M uzerolle et al 



20011) . The apparent NC on top of a BC for i — 15° comes from 
emission near to the disk. It has a base width of the order of the 
typical maximum of 50 K m s of NC in CTTS dBeristain et al 



200ltlBatalha et al.lll996l) . However, the presence of an appar- 
ent narrow component on our model profiles is very dependent 
on the inclination. Such a geometrical effect does not seem to 
be observed. Additionaly, unlike the observed narrow compo- 
nents which show frequently unshifted or redshifted line cen- 



ters dBeristain et alJbOOlLlBatalha et alJI 19961) . this feature has 
a slightly blueshifted center of emission. 

The equivalent width ratios between the different Call 
triplet model line profiles tend to be roughly around the op- 
tically thick regime, usually with the 8542A line slightly en- 
hanced relative to its neighbours. 

It should be kept in mind that the Sobolev approximation 
close to zero velocity does not hold, which is the main reason 
behind the sharp peaks in the line center. In future work, we 
will carry out a more adequate proper escape probability treat- 
ment. Broadening mechanisms and the effects of rotation will 
also be explored. 

3. Model application to Dl Cep observations 

We compare our models to multi epoch observations of DI 
Cep carried out by GFP05 from 1998 to 2001. Further details 
on the observational procedures and instrumentation used are 
given in GFP05. The Call IR triplet lines are blended with the 
Hydrogen Paschen series lines Pa 16-3 (8502.483A), Pa 15- 
3 (8545. 382A) and Pa 13-3 (8665.018A), which are located 
at 157.1 km s _1 , 115.5 km s and 99.6 km s _1 , respectively, 
from the center of the Call lines. There is no indication of 
their presence so we neglect them in our analysis. The pres- 
ence of narrow components in the Call IR triplet line profiles 
is also neglected because the narrow components are not seen. 
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= l.9R e , M* = 



There were no veiling corrections made on the observed lines. 
GFP05 discussed the difficulty of estimating veiling in this 
wavelength region because of the nonexistence or weakness 
of photospheric absorption lines. However, some extrapolation 
can be done from shorter wavelengths leading us to conclude 
that the veiling is low, of the order of ~ 0.05 to ~ 0.2 and in- 
creasing slowly from 8498A to 8662A. In our normalized line 
profiles this implies that the lines can be stronger by 5 to 20 
percent. 

In general, the comparison between our model line pro- 
files and observations should be carried out after removing 
from the latter the underlying photospheric component, i.e, af- 
ter subtracting an appropriately veiled template spectrum. For 
the current DI Cep observations this is not necessary since the 



emission lines are much stronger than the observed absorption 
strength on the template star. 

DI Cep is a classical T Tauri star with spectral classifi- 
cation ranging from F4-K5. The most probabl e spectral typ e 
is between G5-G9 fo r a luminosity class V dlsmailovl 120031) : 
Herbig & Bel] dl988h indicate a G8V s pectral type. Using the 
tables of iKenyon & Hartmannl ( 1995 ) and the in terval G5- 
G9V, the effective temperature is 5590 ± 180K. iKholopov 



(1 19591) determines a distance to DI Cep of 300pc whereas 
Grinin et alj Ql_980) obtain a 200pc distance. We shall consider 
the range 250 + 50pc in o ur analysis. Using photometry from 



Cohen & Schwartz! d 1976b the average < V - I >= 1.41 mag- 
nitu des. Using Kenvon & Hartmannl ( 1 1995b to deredden V - I 



and ICardellietalJd 1989b extinction law (Ry = 3.1) we derive 
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Fig. 8. Contours of x 2 values obtained from model fits to the observations of DI Cep, for varying mass accretion rate and 
maximum temperature of the gas. Darker areas represent better fits. The color/grayscale gradient is not linear, i.e., it its steeper 
near the minimum^- 2 . The white line indicates a confidence level that we set to the upper quartile of the^- 2 probability distribution. 
Each column corresponds to one date of observation. In the last row we plot the intersection of the regions which lie above the 
confidence level on the three previous rows (dots). Extended features in the temperature axis direction are just a consequence of 
model lack of resolution. In the dashed region the continuum importance starts to be non-negligible. 



Ay = 0.88 + 0.42 magnitudes. Using 2MA SS J band photome- 
try an d equation Al from appendix A on lKenvon & Hartmann 
( 1 1995b to derive the bolometric magnitude, yields the luminos- 
ity L* = 3.8 + 1.8L and R+ = 2.0 + 0AR@. One should take 
these values as a mere estimate. Therefore we adopt the value 
R* = 1 .8/? Q and T e fj - 54007^ which results in a stellar mass 
of 1 .42Mt7) ? using the pre-main sequence evolutionary tracks of 
Siessetalld2000l) . 

To determine the inclination i we consider the measurement 



of the period and V sin i found in th e literature. Bouvier et al. 
1986b report Vsnu = 28 km s -1 ; iGameiro & Lagol d 1993b 



found 23 Km s , and in GFP05 Vsin/ is between 19 and 



23 Km s . We adopt Vsin/ = 23.5 ± 4.5Km s . Ismailov 



(12003b found a D I Cep periodici t y from the emission lines 



of 9.24 days and Kol otilov et al.l d2004b 9.2 days from the 



U-V light curve. These spectroscopic measurements result in 
R+ sin i = 4.3 + 0.8/? Q , which largely exceeds R+ derived from 
photometry, and for the adopted stellar luminosity and using 
the same evolutionary tracks, 7?* sin i values imply a mass 
greater than 2M . Given these problems, we adopt the value 
i = 60°, as the models match better the observations. 



We are left with four parameters to set the runs, 
T s hock,Tgas,M acc and the geometry. We have to make as much as 
plausible assumptions as we can to restrict the study to few pa- 
rameters. Therefore, the temperature of the heated photosphere 
in the shock is assumed to be approximately 8000K which 
corresponds to the temperature derived from optical vari abil- 
ity observations b y Gom ez de Castro & Fernandez (1996) and 
Fernandez & Eiroal d 1996b . After an initial exploration of the 



geometry effects on the line profiles, we favoured a large and 
narrow magnetosphere since overall it produces better fits to 
the DI Cep data. The profiles thus produced are narrower, more 
symmetric and exhibit more emission for lower mass accretion 
rates, in better agreement with the data. Independently of the 
values of T s h oc k, T gas and M acc , other geometries are never ca- 
pable of fits as good as those obtained with a large and narrow 
geometry. In particular, the observed variability on DI Cep can 
never be explained in the context of our model by geometry 
changes even when accompanied by changing the remaining 
parameters. Afterwards, we tried to fine tune the mass accre- 
tion rate and the temperature of the gas such that the three IR 
triplet lines could be modeled simultaneously. To accomplish 
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Fig. 9. Best fits of the model to DI Cep observations, in three different epochs for the IR triplet lines. The last two rows in the 
plot show the observed H a and Hp line profiles, for which no models were produced. We have used the gas temperature and mass 
accretion rate corresponding to the regions of Fig.[8]bottom row. 27 July 1999: M„cc = 8.0 x 1O~ 9 M yr _I and T gas = 1 1000K; 3 
August 2001: M acc = 2.5 x 1O- 8 M yr 1 and T gas = 10000 K; 10 August 2001: M acc = 1.6 x 1O~ 8 M yr" 1 and T gas = 12100K. 
Fixed parameters: fl* = 1.8fl , M* = 1.42M , V* = 5400/T, T r = W00K, Geometry = 5.8 - 6fl*, i = 60° 



that we searched our grid of models for the minimum values of 
X 2 , calculated from F b s (observed continuum normalized flux) 
and F,„„d (model continuum normalized flux). This analysis is 
shown in Fig. [8] We plot in white line the confidence level set 
to the upper quartile of the^ 2 probability distribution, on three 
different dates of observation, with each row corresponding to 
fits of the 8498A, 8542 A and 8662 A lines. In the bottom, for 
each line, we plot the intersection of the regions bounded by 
the upper quartile threshold (dots). We also plot the x 2 con- 
tours in the color/grayscale gradient to illustrate the interplay 
between mass accretion rate and temperature. The contours for 
a given line show that the best fits can be produced with M acc 
and Tg as spanning a large range of values. (For example, T gas 
may decrease as long as M acc increases accordingly, and the 
X 2 will remain constant). The darker vertical bands in the x 2 
contours represent areas where an increase on M acc does not 
produce an increase in line emission, keeping x 2 constant for 
the same T gas . The model profiles are saturated which means 
that the lines have completely thermalized. However, as the 
mass accretion rate increases, the assumption of a continuum 
optically thin magnetosphere will eventually become invalid, 
and the optical photospheric spectrum would fade away being 



swamped by an ever increasing continuum emission, putting 
an upper boundary in the x 2 vertical bands. For mass accre- 
tion rates as low as 3 x 1O~ 8 M yr~' in Fig. [8] at a gas tem- 
perature of 12000K, we have to start considering the effects 
from the Paschen continuum emission. For a gas temperature 
of 9000K the magnetosphere starts to have non-negligible con- 
tributions towards the continuum emission for M acc starting at 
about 7.5 x 1O~ 8 M yr . The dashed lines on Fig. [8] mark 
the regions where this magnetospheric Paschen continuum be- 
comes important. We remark that these regions are valid if the 
accretion geometry was approximately axisymmetric as in our 
model. For a non-axisymmetric geometry, with the same T gas 
and same column density, higher mass accretion rates are pos- 
sible since there are fewer accretion columns emitting, and the 
photosphere is also more exposed. 

Using the region of the bottom row in Fig. [8] we made fits 
to the lines on three different observational dates (see Fig. 0, 
illustrating the generally good agreement between model and 
observed line properties. 

Figure [10] shows the equivalent width ratios for all the ob- 
servations of GFP05 (diamonds). Model equivalent width ra- 
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Fig. 4. Source functions for the five transitions considered in 
the Call system, along the streamline that intersects the disk at 
4.8/?*. In the extremities of the magnetosphere the H and K line 
source functions approach each other as well as the infrared 
triplet source functions (see text). We use the same parameters 
as in the low gas temperature case of Fig. 0(c) with the excep- 
tion of the Geometry=4.8-5/?*. Note that a narrower geometry 
has higher densities, so the line source function 8662A is more 
thermalized than in Fig 0(c). 
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Fig. 10. Equi valent width rat i os of t he Call IR triplet lines (di- 
amonds) from lGameiro et aL (2005) observational data sample 
of DI Cep. Equivalent width line ratios for models that repre- 
sent the 27 July 1999, 3 August 2001 and 10 Aug 2001 data 
(Fig.|9]l are shown as asterisks. 

tios obtained from the fits of Fig. [9] are plotted on Fig. [10] as 
asterisks. 

4. Discussion 

The fits to the observations in Fig.|9]are very acceptable, in par- 
ticular given the approximations we are using. Magnetospheric 
infall can therefore be the most important process for the Call 
IR triplet broad component line formation. The Hydrogen pro- 
files show P Cygni character, which arises in a wind. In addi- 
tion, the hydrogen lines are much broader than the Call lines, 



which may also be accounted for by wind emission. Our mag- 
netospheric accretion model is not able to tackle this wind 
emission, and hence the hydrogen line profiles would never be 
reproduced in this context. 

It is possible to model the observed variability by varying 
M aC c and temperature of the gas (Fig. [9]) simultaneously. To fur- 
ther test this conclusion one could do observations/modelling 
of other elements. One could also compare observed absolute 
line fluxes (not available with the dataset discussed here) with 
absolute line fluxes resulting from the models. 

From the results of Fig. [8] and looking at an individual 
line, we can accommodate mass accretion rates ranging from 
5 x 10~ 9 M o yr _1 to 7.5 x 1O~ 8 M yr -1 for a gas temperature 
of 9000K, and a maximum M acc of the order of 3 x 1O~ 8 M 
yr~' at 12000K, for an individual line. Simultaenous good fits 
of the three Call lines for each date of observation enable us to 
restrict the gas temperature and the mass accretion rate to nar- 
row intervals. The variations we find in the gas temperatures 
and mass accretion rates from one date to another lie within 
10000 - 12000K and 7 x 10~ 9 - 4 x 1O- 8 M yr _1 , respectively. 
Taking an extreme upper limit on the accretion luminosity as 
L acc ~ GM*M acc /R t » L, then for the stellar parameters that 
we are assuming, the mass accretion rate must be less than 
M aC c S 1 x lO^ 1 M Q /yr, otherwise L acc would be much higher 
tha n L„. Our values lie comfortably below this limit. 

Gomez de Castro & Fernandez! dl996l) . find that M acc > 



6 x 1O~ 9 M yr~' assuming that all the energy dumped in 
the UV continuum and lines is due to kinetic energy re- 
leased from infall from a corotation radius of 8.3 stellar radii. 



Gomez de Castro & Lamzinl (119991) used the predictions of the 



semiforbidden line ratios SilH/ OIII and Silll/C III for the ac- 
cretion shock (from models by iLamzinl (11998b ). the observed 
semiforbidden line ratios, and the more uncertain values of ex- 
tinction and distance values to DI Cep, to obtain a value of 
mass accretion rate M ar c = 1.9 x 1O~ 7 M yr -1 . They assumed 



the distance of lKholopovl (1 1959b . i.e., 300pc and A v = 0.24. 
Using equat i on 4 o f their work and a distance of 200pc from 
Grinin et al.1 {l980) would bring M acc to about 8.4 x 1O~ 8 M 



yr -1 ; however, using our estimate of Ay one would get 1.3 X 
10~ 7 < Macc % 2.9 x 1O~ 6 M yr -1 , a value which results in an 
acc retion luminosity far too high. 

Johns-Krull eTal] j2000l) find M accl = 4.0 x lO' 1 M Q yr~ l 
and M aC c2 = 7.8 x lO -6 M }>r _1 , using two empirical relations 
between the C/V1549A line excess luminosity with M acc , val- 
ues that would imply, again, an optical continuum veiling much 
higher than observed. We could not determine the influence 
of using our stellar parameters in their work. However, there 
are important sources of uncertainty in their determinations of 
mass accretion rates. In particular, the empirical relations used 
in these determinations are strongly dependent on the extinc- 
tion corrections assumed for the stars upon which the rela- 
tions are derived. In addition, the empirical relations rely on 



Hait i gan et al.l (119951) determinations of mass accretion rates. 



These are based on the boundary layer hypothesis which over- 
estimate the mass accretion rate compared to a magnetospheric 
accretion model. 

Finally, GFP05 find a projected mass accretion rate 
7 M yr -1 using fits of the models of 



M acc cos(i to ) = 2.5 x 10 
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Calvet & Gullbringl dl998l) to their observed wavelength veil- 
ing dependence in the spectral range 4000A - 8000A. Here, 
i t0 is the angle between the accretion tube and the observer di- 
rection. If we use the stellar mass and radius applied in our 
work the projected mass accretion rate would be M acc cos(/ ro ) = 
6.3 x lCr 8 M yr _1 . As the inclination i t0 is not known, the mass 
accretion rate could increase by the factor 1 / cos(/ fo ). 

Given this discussion, it is evident that there is a large 
uncertainty in the mass accretion rate. The values that 
we find are consistent with the l ower limit derived by 
Gomez de Castro & Fernandez! ( 1996 ) and lower than the value 
found in GFP05. The other determinations mentioned are likely 
to be overestimating DI Cep's mass accretion rate. 

The observed Call line peaks in DI Cep follow the relation 

V . > V . > V . (Fig. 15). However, because of the nar- 

8498A 8542A 8662A v & 

rowness of the 8498A line, the equivalent width ratios follow 
the tendency F 8498A ^F 8542A > F 8662A (Fig.[lO). Figure [10] also 
puts in evidence that when the 8498A line becomes stronger 
the others do not have a proportional increase. 

The models reproduce well the line peaks but predict a 
8542A line peak stronger than the other two IR triplet lines 
(Fig. [9). Nevertheless, the more common observed peak pattern 

V . > V .is reproduced. The lines have different intensi- 

JS498A . .8662A , F 

ties mainly due to the different source functions (see figure |4}, 
which follow the relation of the collisional cross sections (<x c ) 
used, i.e., cr c . > cr f . > <x c „. Therefore, the peak ratios 

8542A 8498A 8662A F 

ar e sensitive to the accur acy of the latter, which for example 
in Shine & Linskvl (1 1974b are different from those we use. As 
we are exclusively considering the extended Sobolev approxi- 
mation, and since the line peaks are formed where the approx- 
imation breaks down, i.e., where the material lifts off from the 
disk, our model results around the line center should be taken 
cautiously. 

Figure [10] shows that the model equivalent width ratios are 
of the same order as the observed ones with the exception of 
the more problematic night of 27th of July, 1999. However, 
the variability of the equivalent width ratios does not have a 
complete correspondence in the model. Our variations on M acc 
and T gas alone, do not totally produce the necessary differential 
changes in the profiles to explain their variability completely. 

In some CTTS, specially those with low mass accretion 
rates, the broad component starts to fade a way and a nar- 



row co mponent becomes prominent. From the M uzerolle et al 



(1998) sample we see that the three stars with NC type Call 
profiles have mass accretion rates of about 3.4 x 1O~ 9 M yr _1 , 
which might be interpreted as a lower limit on the mass ac- 
cretion for the disappearance of the BC (M m ,„). In our models, 
depending on the used geometry and temperature of the gas, we 
find in general that M m -,„ should be between 1 x 10~ 9 to 5 x 10~ 9 
M yr _1 , in agreement with the above lower limit. 



An important remark is that lMuzerolle et al.l d 1998b found a 
correlation between the IR triplet line total flux 8542A and the 
mass accretion rate, for mass accretion rates between 1O~ 9 M 
yr _1 and 1O~ 6 M yr _1 . Note that this relation also includes stars 
with NC profiles. Therefore, not only the relation should be ex- 
plained in the broad component regime by our models but also 
the smooth transition between the broad component and the 



narrow component regime should be accounted for. This seems 
to indicate that a multi component model is necessary to pre- 
dict narrow component fluxes coming either from a hot chro- 
mosphere, accretion shock or wind, and the flux es from the BC 



forme d in the magnetospheric accretion flow. iMohantv et al 



(2005) have extended the relation of the Call lines fluxes with 
mass accretion rates to the brown dwarfs. Here M acc can be as 
low as M acc = 1 X 10~ n M yr _1 . Such a multi origin com- 
ponent model for CTTS should also be able to reproduce the 
correlation in this low mass regime. 

5. Conclusions 

The magnetospheric accretion model can explain well the Call 
IR line profiles of DI Cep, as well as the observed variabil- 
ity, using slightly different mass accretion rates and different 
gas temperatures for different epochs. The mass accretion rates 
found for DI Cep that fit the observations are between 5 x 10~ 9 
M yr~' and 2 x 10~ 8 M yr . The peak ratios are not perfectly 
reproduced, although the models can predict the observed ra- 
tio V . > l p .. It is difficult to fully interpret the observed 

8498A 8662A _ J r 

variability using variations of M acc and 7^. Further tests of the 
model and restrictions on the parameter space can be done with 
simultaneous observations/modeling of other elements and by 
trying to reproduce the relation between the mass accretion 
rates and the Call infrared line fluxes. 
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